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Abstract. The Wigner time delay is a measure of the time spent by a particle 
inside the scattering region of an open system. For chaotic systems, the statistics 
of the individual delay times (whose average is the Wigner time delay) are thought 
to be well described by random matrix theory. Here we present a semiclassical 
derivation showing the validity of random matrix results. In order to simplify the 
semiclassical treatment, we express the moments of the delay times in terms of 
correlation functions of scattering matrices at different energies. In the semiclassical 
approximation, the elements of the scattering matrix are given in terms of the 
classical scattering trajectories, requiring one to study correlations between sets 
of such trajectories. We describe the structure of correlated sets of trajectories 
and formulate the rules for their evaluation to the leading order in inverse channel 
number. This allows us to derive a polynomial equation satisfied by the generating 
function of the moments. Along with showing the agreement of our semiclassical 
results with the moments predicted by random matrix theory, we infer that the 
scattering matrix is unitary to all orders in the semiclassical approximation. 



PACS numbers: 03.65.Sq, 05.45.Mt 

1. Introduction 

The dynamics of an open quantum system can be described by its scattering matrix. 
The scattering matrix is defined as the linear operator which transforms an incoming 
wavefunction, expanded in an asymptotic channel eigenmode basis, into the outgoing 
wavefunction. Probability conservation forces the scattering matrix to be unitary. As 
the scattering matrix describes the system, it can be used to investigate the desired 
physical properties of the system. 

A quantity of particular interest is the Wigner time delay [TJ [2] which was first 
derived for the one channel case from a Hermitian operator based on the scattering 
amplitude. The time delay is, as the name suggests, a measure of the extra time a 
particle spends inside the scattering region as a result of being scattered (for some details 
see [3] for example), and was later generalised to multi-channel scattering matrices [4]. 
With M scattering channels, the Hermitian operator admits M eigenvalues which are the 
individual delay times of the system, and the Wigner time delay is simply the average 
of these values. 

For classically chaotic systems with a small opening, the probability for particles to 
remain inside the system decays exponentially with a typical time scale depending on the 
size of the opening. The exponential decay can be seen as a natural consequence of the 
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ergodicity of the classical motion as (for reasonable times) the particle is equally likely to 
hit any part of the boundary of the system leading to a roughly constant small probability 
to leave through the opening each time it hits. The continuous limit of this process is 
the exponential decay, and the time scale associated to the decay is, for chaotic systems, 
exactly the average time delay [5] of the corresponding quantum system. Staying with the 
classical version for now, if we spread particles evenly over the available space and evolve 
them over time they will condense onto a typically fractal pattern arranged around the 
zero-measure set of trapped periodic orbits that never leave the system. From this set, 
and following their stable manifolds very closely we can construct trajectories that start 
outside the system, approach trapped orbits and remain inside the system for arbitrarily 
long times before eventually escaping following the unstable manifolds. In fact, it turns 
out that these trajectories, or rather correlations between them, are responsible in the 
semiclassical limit for recreating the oscillating part of the time delay [6] , which is given 
exactly in terms of the trapped periodic orbits [71 [S] . 

As the oscillating part of the time delay is given in terms of periodic orbits in 
the semiclassical limit (of ft — > 0), correlations between periodic orbits must then be 
responsible for the typical fluctuations of the Wigner time delay. In particular the form 
factor of the time delay can be written in terms of pairs of periodic orbits and the types 
of correlations that contribute were first treated for the spectral form factor for closed 
systems, where Gutzwiller's trace formula [TO] likewise provides a sum over pairs of 
orbits. The semiclassical treatment of this sum started with the diagonal approximation 
[TTj of pairing orbits with themselves (or their time reversal) and used the sum rule of [12] 
to find the semiclassical contribution. This was followed more recently by the treatment 
of correlated pairs of orbits which are almost identical everywhere, but which differ in 
a very small region called an encounter where the orbits behave differently and end up 
reconnecting inside the encounter [13] . All the possible types of orbits with encounters 
were then generated and treated [lU [15] , and this treatment was applied to the time 
delay to obtain a semiclassical expansion for its form factor [16] . 

The expansion for the form factor of the Wigner time delay was shown to agree with 
the result from random matrix theory (RMT) , in line with the idea that most properties 
of quantum chaotic systems are well described by the results of RMT. But RMT can tell 
us at lot more about the expected typical behaviour of the individual delay times, and in 
particular that their distribution should have a compact support [17] in the limit where 
the opening supports a large number of channels. Of course the classical decay gives no 
upper bound on the distribution of the delay times, and so (if we expect to recover RMT 
results) quantum interference as expressed through semiclassical correlations between 
classical trajectories should somehow contrive to limit the maximum delay time. 

Our aim in this article is to identify the sort of correlations that contribute in the 
semiclassical limit and to describe how they lead to the RMT result. We start with an 
introduction to the time delay in section O To study the distribution of the delay times 
we evaluate their moments which can be expressed in terms of correlation functions of 
the scattering matrix. This approach simplifies the semiclassical treatment compared to 
the direct evaluation of the time delay matrix. 

In section 12.21 we briefly review the RMT results before applying the semiclassical 
approximation to the scattering matrix elements to obtain expressions in terms of 
scattering trajectories. For the low order moments, which we treat in section[3J we build 
on the work of the semiclassical treatment of the conductance [TBI [TO] and its second 
moment [20l[2T]. To treat all the moments, we delve into the combinatorial relations that 
first arose in the treatment of the moments of the transmission amplitudes [22]. This 
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work is extended in section 21 where we find an implicit expression for tlie generating 
function of tlie moments and show complete agreement with the RMT distribution of 
delay times. Considering a simpler case in [Appendix A[ we briefly derive the correlation 
functions which recently appeared in the semiclassical treatment of the density of states 
of Andreev billiards [23l. 



2. The time delay matrix 



If we consider a chaotic cavity with one or more open leads that carry M scattering 
channels, the scattering dynamics is encoded in the M x M unitary scattering matrix 
S{E) which relates the incoming and outgoing waves. We are interested in the Wigner 
time delay, which represents the extra time spent in the scattering process compared to 
free motion, and which can be found using the Wigner-Smith matrix [TJ [5J 3] 



i dE 



(1) 



Differentiating the unitarity condition, S'^{E)S{E) = /, with respect to E, we see that 



SHE) 



dS{E) dS^E) 



S{E), 



(2) 



dE dE 

and that the matrix Q is Hermitian with real eigenvalues. The M eigenvalues are the 
individual delay times of the system, and the average value of these times is the Wigner 
time delay 

rw{E) = j^TY[Q]. (3) 
The moments of the eigenvalues are given by 



(4) 



and can be used to recover the complete distribution of the eigenvalues of Q. 

We can also obtain the Wigner time delay from a correlation function of the 
scattering matrix 



sUE-'-^]siE+'-^ 



'2 



(5) 



where it is convenient to specify the energy difference in units of hfi, where fj, is the 
classical escape rate of the system. This correlation function provides a symmetrized 
version of the time delay [23] 



Tw 



1/i de 



2iM 



■Tr 



SHE) 



dS{E) dSHE) 



dE 



dE 



SiE) 



(6) 



which agrees with the definition in (jSj because of the unitarity of the scattering matrix 
as expressed through 

As well as giving us the first moment of the delay times (the Wigner time delay), 
this correlation function can also provide us with the second. Differentiating ([2]) again 
we obtain 



SHE) 



d''S{E) d'^SHE) 



SiE) 



dSHE) dS{E) 



dE^ dE^ ' dE dE 

which can be used to simplify the second derivative of C (e) to 



M 



Tr 



dSHE) dS{E) 



dE dE 



(7) 



(8) 
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Inserting the identity matrix in the form S{E)S'^ {E), we obtain 



de2 



C{e] 



1=0 



M 



Tr [gtQ] = (iAi)2m2, 



(9) 



and the second moment of the delay times. 



2.1. Various approaches to higher moments 

Differentiating the correlation function C(e) further cannot produce higher moments 
since the n-th moment involves the first derivative of S{E) appearing n times while the 
definition of C(e), see (O, involves only two matrices. One way to resolve this difficulty 
is to define a higher correlation function 



C(e,n) = — Tr 
^ ' ' M 



^sUE-'4r\s\E- 



(10) 



where e = (ei,...,e„) is the vector of energy differences. From the rt-th correlation 
function we can obtain all the moments up to the (2n)-th by differentiating. The 
formulae, which are different for odd and even moments, are 



1 



m2k-i 



(i//)(2fc-l) 



(11) 



(iM) 



2k 



(12) 



Though this provides an efficient way to calculate the lower order moments, the 
different energy arguments complicate the semiclassical treatment. On the other hand, 
setting all energy arguments to e and differentiating n times leads to appearance of 
additional unwanted terms. To deal with such terms we define the correlation function 



D(e,n) = — Tr 
^ ' ' M 



sUe~'^\s(e+'^\-i 



[/(^)-/(o)r 



e=0 



[/'(O)]", 



Because of the relation 
1 d" 
n! de" 

we can see that the moments of the time delay matrix are now given directly by 

1 d" 



(i^)"n! de' 



:D{e,n) 



Expanding the n-th power in (|13l) we can also obtain the moments as 



1 



(i^)"n! de" ^ 



(-1)' 



i-k n 



in terms of the correlation functions 



sUe 



ehfi 



S E 



ehjjL 



(13) 



(14) 



(15) 



(16) 



(17) 



The correlation functions C(e,n) have applications outside the scope of this article 
(see [33]) and thus their semiclassical evaluation is a question of stand-alone interest. 
However, we found that performing the summation in (|16p is a difficult task. Thus 
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from now on we will work with D(t, n) directly. However, the technique we develop for 
evaluating D{e,n) is suitable for the simpler function C{e,n) as well. We exploit it by 
deriving an equation for the generating function of C{e,n) in [Appendix A[ 

Finally, we note that one can use the semiclassical approximation directly in 
the definition of the Wigner-Smith matrix, equation ([1]). The approach we took, 
however, allows us to build on previous semiclassical work performed in open systems, 
in particular on the average conductance of a chaotic ballistic device and its moments 
[m EHl EU [22] ■ We will obtain, as in the case of the conductance and shot noise, simple 
diagrammatic rules for the semiclassical contributions of correlated trajectories. But 
before approaching this task we will quickly review the RMT results for the delay times. 



2.2. Random matrix predictions 

The random matrix result for the probability distribution of the delay times is |17| 

Pir)^:^Vir+-r) {t~t^), r± = (18) 

from which we can calculate the moments 



J T^pir) dT^^jJ T"-'Vir+-r) {r - r^) dr. 



(19) 



Substituting t + 3 — fir this can be written more obviously in terms of the moments of 
the semicircle distribution 

m^^^ r {t + Sr-' ^'^^l;'" dt, (20) 



M" J-R T^R^ 

with = 8. The odd moments of the semicircle distribution are 0, while the even 
moments can be written as 



'^^^ = I J d^=(-) c„ = 2"c„, (21) 



where Cn are the Catalan numbers defined as 
1 f 2n\ (2n)! 



(22) 

n + 1 V J (n+ l)!n!' ^ ' 

By expanding the term {t + 3)"^^ in (j^U)) in powers of t, we can therefore express the 
moments (beyond the first) by the following sums 

^ k=0 ^ ^ 

--+3^-^E(';I>-3^^"-'=^2'=^- (24) 
^ fc=0 ^ ^ 

The sequence of integers obtained from combining the Catalan numbers as above, and 
hence appearing in the moments, can then be shown |25| to be the so-called Schroder 
numbers (Sloane's A006318 or 2x A001003) [26]. More importantly for our purposes, the 
Schroder numbers can be derived from a generating function '26' and so the moments 
of the delay times can also be generated from 
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3. Calculation of the low correlation functions 

Since the matrix S{E) is unitary we can rewrite the correlation function D{e,n) as 

' ' ■{E-'-!l^)s[E+f)-SHE)SiE)f (26) 

1 E. (27) 



where the summation in the second hue is over choices of n incoming channels and 
the ± subscripts represent the different energy arguments of the scattering matrices. 
To evaluate D(e, n) semiclassically, we will use the semiclassical approximation for the 
scattering matrix elements, which is given in terms of open trajectories |27[ I28[ 118] 

1 



SUE)-^ Ace^^^. (28) 



The Th appearing in the prefactor is the Heisenberg time, and it is simply related 
to the classical escape rate hy fi — M/Tn. In the sum which is over all classical 
trajectories C, that start in channel a and end in channel h (where the channels fix the 
absolute value of the angles at which the trajectories enter and leave the cavity) 5*^ is 
the action of the trajectory C, and its stability amplitude [35] including the phase due 
to the number of conjugate points along the trajectory. As our starting point, we will 
substitute approximation ((28|) into ((26)) . expand the action up to first order in energy 
Sq{E + 5E) K. Sc_{E) + T^{E)5E, where T,^ is the time the trajectory ( spends inside the 
system, and ignore any change in the slowly varying prefactor . 

Below we perform this calculation for n = 1, 2 and 3 before formulating the general 
counting rules and performing the evaluation of D(e,n) for general n. 

3.1. Calculating D(e,l) 

After semiclassical approximation (|28|) . the correlation function D{e, 1) becomes 

which is a sum over trajectory pairs C' both of which start and end in the same channels 
(ii and oi respectively), followed by a sum over all the possible channels. We note that 
the difference of actions S(; — S(;' is divided by ft ^ 1 and therefore the resulting phase 
oscillates wildly unless the action difference is of order h. The semiclassical expansion 
is based on identifying couples (or, more generally, families) of orbits that have small 
action differences. 

To leading order in inverse channel number 1 /M the first correlation function can be 
calculated as in [B] using the diagonal approximation This approximation restricts 
the sum to trajectories C and C that are identical 

^"''(''^) = lik.^ ^ lAcP(e-^^^-l). (30) 

il.oi C(il— i-oi) 



Moments of the Wigner delay times 



7 



The sum in (j30l) can be performed by using a sum rule for open trajectories |18] which 
turns it into an integral over the trajectory time T 

J2 \Ac\'-.-^ TdTe-'^^.... (31) 

C(il->oi) 

The exponential term in (I3ip represents the average probability that a trajectory remains 
in the system for the time T, while the sum over channels, where we can pick both 
ii and oi from the M possible channels, simply gives a factor of Af^. The diagonal 
approximation thus gives 

D^'-s(^e, 1)-^ p AT _ ^-^t^ (33) 

1 = T^, (33) 



1-ie 1-ie 

where we have used the fact that /i = M/Tn to simplify. Substituting into (jlSp we 
obtain 

rw = (34) 
M 

which just states that the average delay time is the inverse of the classical escape rate 
as we might expect. 



3.2. Calculating D{e,2) 

To obtain the next moment we move to the correlation function Z?(e, 2). For this we 
will need to treat not only the diagonal pairs, but also correlated trajectories that have 
encounters. The calculation follows the calculation of the shot noise power [50] which 
was first performed to leading order in inverse channel number for quantum graphs [29] . 
The semiclassical treatment of the shot noise builds on work on the conductance [HI [19] , 
which itself is built on the work on spectral statistics [TS] [M] [15] . 

First we slightly modify (|27|) to explicitly include an indication of the unitarity of 

S 

D{e,2) = ^ y^(slS+~5,,,,,S^s) (sls+ - S,,,,,S^ s) . (35) 

While the semiclassical approximation (|28|) preserves the unitarity of S, and thus the 
Kronecker deltas are not necessary, their inclusion greatly facilitates the derivation. 

We now write down the semiclassical expression for the correlation function D{e, 2) 
in terms of open trajectories 

^(^'2)«^5: Y: E Ac^c'A,A|,e^(^c-5,+5,-S,) 

^ il,i2 C(*l-i-oi) 4(12^02) 

C'fe^Oi)C'(4i->02) 

X (e^(^<+^c') - S,,,^ [e'-^W,) _ ^^^^^^^ (36) 

For diagonal terms we can either pair ( — C and ^ = or ^ = ^' and ^ = C'- For the 
first case, the start channels ii and 12 must coincide thus triggering the Kronecker delta. 
For the second case we have oi ~ 02- Either case leads to a (leading order) channel 
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Figure 1. An example of two trajectories with a single encounter and two partner 
trajectories. 

factor of , because the number of channels we can choose in the outer summation is 
3. Putting these two possibihties into p6p we can write the diagonal contribution as 

2) = ^ ^ \A,\Mi? (e'^^^^ - l) (e'^^^^ - l) 

+ ^El^clVdV^^^^e'^'^^^ (37) 

Note that in the second case it can also happen that ii — Z2, but this case gives a 
lower order channel factor of and is therefore neglected. Using the open sum rule 
from (|3T]) we get 

i^d-g( 2) = -^^^ + L^. (38) 

(1-i.)^ (1-ie)^ 

However, as we know from the calculation of the shot noise [291 ED], the diagonal 
terms are not the only ones that contribute to leading order in inverse channel number. If 
the trajectories C, and ^ come very close to each other in an encounter, as in Figure[TJ then 
the partner trajectories can cross over inside the encounter leading to a quadruplet of 
trajectories with a small action difference. Such a quadruplet can then give a contribution 
in the semiclassical limit. While such an encounter makes the contribution higher order 
in inverse channel number, the sum over channels now contributes the factor of Af^. As 
a result the quadruplet contributes at the same order as the diagonal terms. 

To calculate the contribution, we simply put the additional energy arguments into 
the calculation of the shot noise power. The contribution can be separated into a product 
over the links and the encounters |20j and written as 

1 m3 r°° r°° 

Di^'){e,2)^ ^ dti e-''(i-'^)*W dta e-^^^-'^)*^ 

-'O Jo 

oo 

di. e-''(i-'^)*^ / dt4 e-^(i-'^)*^ 



3e / clt4e 

Jo 

p-Ai(l-2ie)to„c(s,M)p7jSii 

dsdu — ^ , (39) 

S2ienc(s, U) 
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where (2^) refers to the structure of the diagram: one encounter with two (unprimed) 
trajectories meeting (a "2-encounter" ) . In the final integral, s and u are the separations 
along the stable and unstable manifolds of the two original stretches inside the encounter 
and f2 is the volume of the available phase space, while tonc(s, u) is the duration of the 
encounter. An important point is that although the encounter involves two trajectory 
stretches, as they are close to each other they will either remain inside the system 
or escape together and their average survival probability is given by the time of just 
a single stretch. The presence of the encounter actually slightly enhances the survival 
probability of the whole trajectory quadruplet, and this tiny classical effect has important 
semiclassical implications. Elsewhere in (p9l) . the ti are the durations of the link stretches 
as depicted in Figure [TJ Performing the integrals following [20] , we obtain the result 

= (40) 

The structure of the answer is very simple: each Z-encounter contributes a factor of 
— (1 — iZe), while each link stretch gives the factor (1 — ie)"^. These diagrammatic rules, 
which first arose for the conductance ^19 , massively simplify calculating the semiclassical 
contributions, and are the reason we consider correlation functions rather than the time 
delay directly. 

We can combine the two leading order results from (p8| and (|40|) to obtain the 
second moment 

2 , , 

ma = (41 

M 

In conclusion of this subsection we mention that the configurations described as 
diagonal contributions above can be obtained from the diagram in Figure [T] by setting 
ti = t^ = Q for the first case and t2 — — Q for the second. We will refer to this 
reduction as moving or sliding an encounter into the lead. As we have seen above, 
moving an encounter into the input lead can lead to a contribution that is different from 
the encounter in the output lead. 



3.3. Calculating D{€, 3) 

Before we proceed to calculate D{e,3), we will briefly look at how we can form the 
diagrams that contribute. As we have seen, the n-th correlator is expressed as a sum 
over 2n trajectories 

{Cj(jj+l^Oj)} 

(42) 

Taking the trace of the product of matrices means that we identify in+i — h- Therefore 
the trajectories complete a cycle, if we consider moving forward along the unprimed 
trajectories and back along the primed ones. The resulting structure for n = 3 is shown 
in Figure [2^ and, as we have also seen, in (|42|) we add the actions of the unprimed 
trajectories and subtract the actions of the primed ones, so the resulting phase oscillates 
wildly unless the total action difference is of the order of h. 

To obtain such a small action difference we can collapse all the trajectories onto 
each other, as in Figure ^p, creating encounters of the type we saw in Figure [T] It 
turns out that to obtain all contributions from this type of direct collapse we need to 
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Figure 2. The original trajectory structure for D{t, 3) in a) can be collapsed down 
to (tliree copies) of tlie structure in b), which has two 2-encounters, to ensure a small 
action difference. Sliding the encounters together creates the single 3-encounter in 
c), while moving encounters into the leads generates all the further possibilities. 
Trajectories C,j are indicated by solid lines and trajectories by dashed lines. 



cyclically permute the labels of the trajectories, resulting in three copies of the diagram in 
Figure Od. Alongside this direct collapse we can imagine sliding the encounters together 
to create a single diagram with a single 3-encounter, as shown in Figure [2j;. Further 
possibilities then arise from sliding the encounters into the leads, giving the remaining 
diagrams depicted in Figure [51 

Before we write down the total contribution, we consider two instructive examples. 
First we evaluate the contribution from the diagram in Figure [5^. We assume the 
trajectories indicated by solid lines are numbered top to bottom as d, C2 and Cs- Then 
we have ii = 12 since the first encounter is in the incoming lead. This activates one of 
the Kronecker deltas in the expression for D(e, 3) (see or the similar expression in 
([55]) V The terms due to trajectories C2 and Cs, which have a non-degenerate encounter, 
are the same as in ([5^ . giving in total 

-,(2") Q^ _ n Arr. ^-t^T., (Jet^Tc, _iA TT TT n A^rP „-^(i-'^)'^c. 



J,r4 /"CO pco 

-'0 j=2,3p=l,2-^0 



e 



g--Ai(l-2ie)tc„c(s,ti)g-^sii 



where refers to the duration of p-th part of trajectory Q and the power of M came 
from the 5 choices of the remaining channels. Evaluating the integrals we get 

If on the other hand, we consider the contribution of the diagram in Figure [5^, we 
notice that 11=^2= is, thus activating three Kronecker deltas. The contribution of this 
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diagram is thus 

^^..(^.3) = f ri/°"dTc, e--<. (e--<. - l) ^ (45) 

Comparing the results in (03]) and P5)) we can surmise that the power of ie 

that results from sliding an encounter into the incoming lead is equal to the number of 
direct stretches from the encounter to the outgoing lead, i.e. the stretches that do not 
participate in any other encounters. This observation can be mathematically verified 
by using the relationship between diagrams and factorisations of the cyclic permutation 

EH]. 

To summarise, we can write down the contribution of each diagram by simply 
looking at its links and encounters and assigning 

• a factor of (1 — ie)^^ to each stretch, 

• a factor of —(1 — ile) to each non-degenerate ^-encounter, 

• a factor of (ie)'* to each encounter happening in the incoming lead and having s 
direct stretches to the outgoing lead. 

• a factor of 1 to each encounter happening in the outgoing lead. 
Altogether, the leading contribution to D{e. 3) is thus 

'{l-2\ef -(l-2ie) -ie(l - 2ie) ie 



£)(e,3) = 3 



(l-ie)7 (l-ie)5 (l-ie)5 (1 - ie) 
-(l-3ie) 1 (ie)3 



(46) 



(l-ie)6 (l-ie)3 (l-ie)3' 

where the factor of three counts the three different ways to label the diagram of Figure[2]3 
and its descendants. 

The third moment can thus be calculated to be 

^3 = (47) 

and we note that if we put in three different energy arguments in line with ([T0| , we can 
use the same diagrams to get the next three moments using ([TT|) and as 

22 90 394 

"15 = -T, TO6 = — (T- (48) 
/i"" /i'' ^jP 



4. All moments 



Now that we know the rules which govern the contributions of individual diagrams, we 
can look at generating all the diagrams and their contributions recursively. To leading 
order in inverse channel number, the contributing diagrams of the type 2"^ 3"^ • • • are in 
a bijective correspondence with rooted plane trees [2^ that have vi vertices of degree 
21 and all other vertices of degree one (called "leaves" ) , see Figure [H We denote by 
V = V2 + + . . . the total number of vertices of degree higher than 1. These vertices 
correspond to encounters in the diagram. The total number of leaves can easily be seen 
to be 2(L — + 1), where L = 2v2 + 3v3 + . . .. Starting from the root the leaves are 
labelled ii, oi, 12, 02, . . . , in, o„, where n — L — V + 1 is. the order of the correlation 
function D{e,n). The i-labelled leaves correspond to incoming trajectories starting in 




the lead and the o-labelled leaves correspond to the outgoing trajectories exiting into 
the lead. 

As we have seen, some encounters can touch the lead, but this can only happen if 
an ^-encounter has I vertices with labels i connected to it («-touch) or I vertices with 
labels o connected to it (o-touch). For example, in Figure [31 the vertex a can i-touch, 
the vertex c can o-touch and the vertex b can do neither. 

When an encounter can touch the lead, the total answer we seek is the sum of the 
contributions from when it does and when it does not do it. Equivalently, we can take 
the multiplicative factor of an encounter to be the sum of all factors it can produce. To 
illustrate this point, we revisit the calculation of D{e, 3) and rewrite equation in the 
form 

n( ^^ 3 ^ -(l-2ie) , ^ -(l-2i6) 

1 /-(l-3ie) 



This is the sum of contributions of three (2^)-trees and one (3^)-tree, see FigurelH The 
structure of a contribution is as follows: the prefactor is (1 — ie)~", where n is half 
the number of leaves (the order of the correlation function). Then follow the factors 
corresponding to the vertices of the diagram. Each vertex of degree 21 (corresponding 
to an Z-encounter) gives a multiplicative factor of 

71— )r' (50) 

modified by an additional +1 if the vertex can o-touch and by +(ie)'* if the vertex can 
i-touch, where s is the number of o-leaves attached to the vertex. In the case when there 
is only one vertex in the diagram, as in the (3^)-tree in Figure |4l it can both i-touch and 
o-touch but not at the same time. To provide a further example, the overall contribution 
of the tree in Figure [3] is 

To count all possible trees while keeping track of the structure of their vertices 
we introduce the generating function F{x, z'o, zI,t). The roles of the variables are as 
follows: 
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Figure 4. The trees contributing to the third order correlation function 3). 
The first three trees correspond to relabellings of Figure [2J3 and the last one to 
Figure [2};. 




Oi 12 02 h 03 k O4 is O5 

Figure 5. A separation of a tree into the top vertex and subtrees. Note that in the 
middle subtrees the roles of variables (zi,a) and {ii,T) are reversed. 



• the power of xi enumerates the number of non-degenerate ^-encounters 

• the power of Zo^i enumerates the number of Z-cncountcrs that o-touch the lead 

• the power of 1 + cr is the total number of i-labelled leaves adjacent to the encounters 
that o-touch the lead 

• the power of Zj / enumerates the number of ^-encounters that j-touch the lead 

• the power of H- r is the total number of o-labelled leaves adjacent to the encounters 
that i-touch the lead. 

For example, the tree in Figure [3] gives rise to four contributions to the generating 
function, corresponding to the four possibilities of the vertices a and c touching the lead 
or not 

X2X3X2 + ^1,2(1 + t)x3X2 + X2X3Zo.2{'^ + 0-) + Zi,2(l + t)x3Zo,2{.'^ + Ct) 

- (x2+ Zi^2{l+T)^X3(x2 + Zo,2{l + (j)y (52) 

Our aim then is to set 

^i^~rr — ^^T^' Zo,i^Zi,i^l, cr = 0, r = -H-ze, (53) 

in line with the semiclassical contributions described above, and finally include a change 
of variables to provide the correct prefactor of (1 — ie)~". 

While the generating function F is our aim, we will first deal with an auxiliary 
function / = f{x,Zo,(J,Zi,T) which is defined exactly as F except for not allowing the 
top vertex to i-touch the lead and not counting the root as a leaf when the top vertex 
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o-touches the lead. These restrictions on the function / makes it possible to find a 
recursive equation for it. The value corresponding to an empty tree will be set to 1, 
i.e. /(O) = 1. The correlation function calculated with the above restrictions will be 
denoted by D{e,n). 

To write a recursion for / we separate a tree into the top vertex of degree 21 and 
2/ — 1 subtrees, see Figure [5] If the top vertex is non-degenerate, its contribution is 
xif^ , where the function / correspond to the even-numbered subtrees in which the 
positions of z's and o's are reversed. Thus the roles of all the variables corresponding to 
leaves of one type are switched, i.e. / = /(x, i^, r, i^, a). 

To account for the possibility of the top vertex touching the lead, we recall that from 
the definition of /, it is only allowed to o-touch. In this case all odd- numbered subtrees 
must be empty and the contribution of each even- numbered subtree is f + t. Note that if 
the subtree is empty, the contribution evaluates to the correct value of /(O) + t — 1-1- r. 
Putting this together, we have 



f = l + }Axiff-'+ZoAf + r) 



and correspondingly 

/ = 



1 



OO 

oo 
1=2 



Xlff- 



(54) 



(55) 



To obtain the function F we need to take into account the special role of the top vertex. 
It can both i-touch and o-touch although not at the same time. Additionally, it is always 
adjacent to an i-labelled leaf and thus always contributes to the power of 1 -I- t. The 
final generating function F satisfies 



F = 1 



OO 

E [^iff'^' + ^oAf + ry-'{l + r) + zuif + 'jY 



1^2 



oo 



1=2 



1=1 



[rZoAf + ry-^ + ^z,Af + ^)'"' + fz^,l{f + ^)'~' 



(56) 



where, to get to the last line we defined Zo,i — Zi^i = 1. We note that the value F{0) 
cannot be defined from recursive considerations and needs to be chosen to provide the 
correct answer for D(e, 1), which turns out to be ie = 1 + ct + t. Another important 
observation is that our choice of the leaf ii as the root (see Figure [3] for example) is 
arbitrary. In particular an o-leaf could be chosen and the answer for D{e,n) should 
not depend on the choice. Thus the function F should be symmetric with respect to 
swapping the variables Zq with Zi, r with a and / with /. This is not apparent from 
([55)1 but will be checked (and used!) at a later stage. 
Now we can make the substitutions 

- ' (57) 



(1-ie)' 



(1 
0, 



le. 



(58) 
(59) 

(60) 
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The substitutions give the contribution of each tree diagram as in (I53p . but we 
have included the powers of r to keep track of which order correlation function they 
contribute to. Indeed, the power of f corresponding to a 2'"'^'S"^ ■ ■ ■ tree would be 
W2(2 — 1) + U3(3 — 1) + . . . — L ~ V — n — 1. Substitution ((5^)) therefore gives a 
prefactor of (1 — ie)^^" and so to get the additional factor of (1 — ie)~^ we need for the 
correct prefactor, we introduce g = (1 — ie)~^/. The function g is a generating function 
of the "restricted" coefficients D{e,n), i.e. 

oo 

g = J2r"-'D{e,n). (61) 

n=l 

Performing all the changes of variable, apart from (I60p for now, from (|54p we arrive at 



9 



(l-ie) 



i-i 



(62) 



g{l ^ie)^l-J2 r'-'g'g'-'il He) + J] r 

1=2 1=2 

and a similar equation for g. The sums can be performed easily, especially when we 
notice that the first two terms correspond to the I — 1 terms of the two sums, leading to 

9 ieff , 1 



1 - rgg (1 - rggf 



1 



rg - 



(l-i£) 

Taking the numerator of the above equation, we arrive at 

1 + rgr 



(1 - ie - rgg) 
and similarly for g 

(1 - ie - rgg) 



1 



rga 



+ iieyrgg = 0, 



+ {iefrgg = 0. 



(63) 



(64) 



(65) 



We have written the equations in a form that highlights the symmetric terms involving 
gg. Taking the difference between and we obtain 

g{\ — ie — rr) = g(l — ie — ra). (66) 

This can now be substituted back into ([M)) or (p5)) to give an implicit equation for g or 
g, although we only state the result after the simplifying substitution from ([50]) 



9 



(1-ie) 



(1 + 



rg 



1 



(1-ie) (1-ie) 



(67) 



(1-ie) 

Putting all substitutions apart from (1501) into (|56l) . and defining G = r(l — ie)^^F, 
we obtain 

G - - (68) 



1 



rg 



(1 — rg)(\ — ie) — rcr (1 — rg)(l — ie) 



(1-ie) 

Though it is clear that the last two terms together are symmetric, this still needs to be 
checked for the first term. We can rewrite it as 

5(1 -ie) _ 55(1 



9 



1 



rcr 



rg{l 



.9(1 



le 



ra) - rgg{l - ie) ' 



(69) 



so that we can see the symmetry follows from (|66p . 

Having verified the symmetry of the generating function G we can put in the 
symmetry-breaking substitution (j60p . or rather the simpler and equivalent r = and 
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(T = ie — 1 . We finally get the generating function of the required semiclassical correlation 
functions 

G-T^%7^' G-Y.^-D{e^n). (70) 

' n>l 

From (jlSp . we can see that the n-th moment is simply (i/z)~" times the coefficient 
of (er)" in the expansion of G. Setting a = ie to cancel the factors of i, we just need to 
extract the coefficients of (ar)" from G. To do so we let y — r{g — 1) and obtain 

G = -^, y = j^- (71) 

1 — y 1 + G 

Substituting g = l + y/r, r = s/a and setting a = in equation (|67|) . we get the equation 
s~ys~y + 2y^ = 0, (72) 
for y and, through relation ([7T|) . the equation 

G^ + {s - 1)G + s = 0, (73) 

for the generating function of the moments. Solving this, and taking the solution which 
gives the correct value of G = when r = 0, we get exactly the generating function 
(|25p from RMT. As the semiclassically calculated moments match the RMT ones to all 
orders, we then recover the full distribution of the delay times p8)) from |17j . 



5. Conclusions and Outlook 



For open chaotic systems we presented a semiclassical derivation of all moments of the 
individual delay times, up to leading order in inverse channel number. The derivation 
essentially relies on the twin properties of hyperbolicity and ergodicity of the chaotic 
dynamics: the latter for the possible return to nearby points to create encounters and the 
former to allow the reconnections inside those encounters. This leads to a result which is 
in agreement with the RMT prediction for the delay time distribution, a semicircle type 
law. Notably this implies an upper bound on the longest delay time. In the derivation 
we relied heavily on the previous work on semiclassical expansions, in particular jl9[|21| . 
which contain an implicit assumption of instant equilibration in the underlying classical 
dynamics. The influence of slower equilibration can be explored by treating the effect 
of the Ehrenfest time on the semiclassical contributions [23l [31] . When the Ehrenfest 
time becomes much larger than the typical time trajectories spend inside the system one 
recovers the classical exponential distribution of delay times. 

It is also interesting to explore the effect of moving from chaotic to mixed phase 
space (with regular islands) on the distribution of the delay times and this upper bound. 
The interference effects which lead to the current result rely on the chaotic dynamics 
and should be suppressed if the chaotic part is reduced. However, additional effects such 
as periodic orbit bifurcations [32] can also be fairly strong. The moments of the delay 
times and in particular their upper bound could therefore be very sensitive measures to 
explore the dynamics inside quantum dots, and possibly used to measure the relative 
weights of the chaotic and regular parts of phase space. 

By considering scattering matrix correlation functions, which were then treated 
semiclassically, we were able to derive equations that implicitly define the generating 
functions of the moments. It is worth noting here that such correlation functions are 
useful for investigating other questions, such as the density of states of chaotic Andreev 
billiards. 
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Restricting our attention to the leading order in inverse channel number resulted 
in limiting the contributing diagrams to trees only. This, in turn, allowed standard 
recursive tools to be used. Looking beyond the leading order, the semiclassical diagrams 
have more complicated structures so that such tools can no longer be directly applied. 
Thus, in our view, the remaining challenge of obtaining a semiclassical expansion of all 
moments to all orders is a task of significant technical difficulty. However, it is a task of 
particular interest, not only because the subleading orders should be influenced by the 
symmetries of the system, but also because the effects of a finite number of channels is 
of much experimental relevance. A solution to this problem would have to involve new 
combinatorial tools and we hope that a clear and general algebraic structure will emerge 
as a result of research in this direction. 
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Appendix A. Correlation coefficients C(e, n) 



Here we briefly outline the results for the correlation functions C{e,n) defined in (ITT)) . 
As these do not involve subtracting the identity matrix from each bracket (as we did in 
([13])) we no longer need to change the contribution of certain diagonal pairs or subtract 
anything when the start channels coincide. We can therefore simply set both a ~ t = 
in the treatment of section 2] and generate C(e, n) instead. This simplification means 
that g — g because of (j66|) . so that g is now given implicitly by 



1 



1 



rg 



;i-ie) 



1 



(l-i6) (l-i6) 



(A.l) 



where in fact the only difference is that the r on the left hand side of ((67|) has disappeared 
We can generate the first few terms in the expansion <? = 



1 



30 = 



31 



92 



g„ as 
{2e- - 1) 



(A.2) 



(1-ie)"' (1-ie)' 
Adding the contribution from the top node in F, we obtain the full generating function 
of C(e,n), which is given by 

rg 



G{e,r) 



rg + r^g"^ + 



1 — rg 

The expansion G{e,r) — Xr"C(e,n) then leads to 

1 



C{e, 1) - go = 
C(e,2)= [gi 



(l-i6)' 

+ .90^] = 



C7(e,3) - [52 + 2gigo+3o'] - 



1 - 2ie - 2e^ 

1 - 4ie - 9e2 



(1 



(A.3) 

(A.4) 
(A.5) 
(A.6) 



which are exactly the correlation functions we can obtain by considering the diagrams 
explicitly as in section |3l Combining these results in line with expanding (|13p . we 
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recover the first three functions D{e,n) calculated in section [3] We could also continue 
to generate terms to obtain the moments via ((T5)) . 

However, the function G(e, r) contains more information than just the moments, for 
example by setting e = we see that g = I- Hence 

G(0,r) = = r + r2 + ..., (A.7) 
1 — r 

which shows that 

Tr[SHE)S{E)]" ^ M, (A.8) 

and that the unitarity of the scattering matrix holds semiclassically for all powers n 
to leading order in inverse channel number. As another example, this function G and 
the correlation functions it generates appear in the density of states of chaotic Andreev 
bilhards (23) . 
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